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ABSTRACT 

The sensitivity of pulsar timing arrays to gravitational waves is, at some level, limited 
by timing noise. Red timing noise - the stochastic wandering of pulse arrival times 
with a red spectrum - is prevalent in slow-spinning pulsars and has been identified 
in many millisecond pulsars. Phenomenological models of timing noise, such as from 
superfluid turbulence, suggest that the timing noise spectrum plateaus below some 
critical frequency, / c , potentially aiding the hunt for gravitational waves. We examine 
this effect for individual pulsars by calculating minimum observation times, T m i n (/ C ), 
over which the gravitational wave signal becomes larger than the timing noise plateau. 
We do this in two ways: 1) in a model-independent manner, and 2) by using the 
superfluid turbulence model for timing noise as an example to illustrate how neutron 
star parameters can be constrained. We show that the superfluid turbulence model 
can reproduce the data qualitatively from a number of pulsars observed as part of 
the Parkes Pulsar Timing Array. We further show how a value of f c , derived either 
through observations or theory, can be related to T m j n . This provides a diagnostic 
whereby the usefulness of timing array pulsars for gravitational-wave detection can be 
quantified. 
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1 INTRODUCTION 


Pulsar Timing Arrays (PTAs; e.g., Manchester et aI7||2013 


|Kramer fc Champion|2013[ McLaughlin|2013 1 seek to detect 

nanohertz gravitational waves from cosmological and extra- 
galactic sources by looking for correlations between contem¬ 
poraneously measured pulse arrival times from multiple ra¬ 


dio pulsars (Hellings & Downs 19831. The sensitivity of a 


PTA is limited by pulsar timing noise, i.e., stochastic wan¬ 
dering of pulse arrival times. External noise sources include 
interstellar plasma turbulence, jitter noise and errors in ter¬ 
restrial time standards; see Cordes (2013) for a description 


of all dominant noise sources and an estimate of their mag¬ 
nitudes. Intrinsic noise sources have been attributed to mi¬ 
croglitches (Cordes & Downs 1985 D’Alessandro et al.|1995 


et al. 2006 


Melatos et al.|2008l, post-glitch recovery (Johnston & Gal 


loway||1999 1, magnetospheric state switching (e.g., jKramer 


_ Lyne et aL|2010 l, fluctuations in the spin-down 

torque (Cheng 1987 b|a Urama et al. 2006), variable cou¬ 
pling between the crust and core or pinned and corotating 


regions (Alpar et al.|1986 Jones|1990 1, asteroid belts (Shan 


E-mail: paul.lasky@monash.edu 


non et al. |2013[) and superfluid t urbulence ( |Greenstein|1970'} 
Link 2012 Melatos & Link 20141. Analyses of long-term mil¬ 


lisecond pulsar timing data indicate that timing noise power 
spectra are typically white above some frequency and red 


below it (Kaspi et al. 

1994 Shannon & Cordes 2010 van 

Haasteren et al. 2011 

Shannon et al. :2013). 


Red timing noise power spectra cannot extend to ar- 
bitragerily low frequencies, as the infinite integrated noise- 
power implies divergent phase residuals and hence (if phase 
residuals arise from torque fluctuations) unphysical pulsar 
angular velocities. One therefore expects the spectrum to 
plateau, or even become blue, below some turn-over fre¬ 
quency f c - A number of physical models naturally pre¬ 
dict low-frequency plateaus, including superfluid turbu¬ 


lence (Melatos & Link 20141 and asteroid belts (Shannon 
et al. 2013:. We discuss the former in detail below. A low- 
frequency plateau enhances prospects for the detection of 
a stochastic gravitational wave background. As the gravita¬ 
tional wave spectrum is a steep power law for most cosmolog¬ 
ical sources (e.g., Maggiore 2000, Phinney 2001; Grishchuk 
20051, it rises above the plateau below some frequency as 
long as it too does not have a low-frequency cut-off (e.g., 
Sesana|2013a| Ravi et al.|2014 and discussion below). 
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2 Lasky, Melatos, Ravi & Hobbs 


In this article, we quantify how a low-frequency tim¬ 
ing noise plateau affects the direct detection of gravitational 
waves with PTAs. Specifically, we calculate the minimum ob¬ 
servation time for any individual pulsar to become sensitive 
to gravitational wave stochastic backgrounds from binary 
supermassive black holes (SMBHs) and cosmic strings. We 
note this minimum observation time is only an indicative 
quantity for determining when a gravitational wave signal 
will dominate the timing residuals for an individual pul¬ 
sar; it does not account for algorithms that correlate noise 
properties between pulsars, a point we discuss in more detail 
throughout. We do this in two ways, firstly by parametrising 
the timing noise in a model independent way, and secondly 
by applying the superfluid turbulence model of |Melatos fe| 
Link (20141. In the first approach, we express this minimum 


time in terms of three pulsar observables; the amplitude and 
spectral index of the timing noise power spectral density and 
the turn-over frequency. The second approach is included as 
an example of how to relate PTA observables to neutron star 
internal properties in the context of one particular physical 
model with only two free parameters. It does not imply any 
theoretical preference for the superfluid turbulence model, 
and will be extended to other physical models in the future. 

The paper is set out as follows. In section [2] we define 
a phenomenological model for timing noise, and review pre¬ 
dictions for the power spectral density of the phase residu¬ 
als induced by gravitational waves from SMBHs and cosmic 
strings. In section[3]we calculate minimum observation times 
for hypothetical pulsars as a function of their timing noise 
spectral index, normalisation, and turn-over frequency. In 
section [4] we apply the superfluid turbulence model to data 
and extract ‘by-eye’ parameter estimates for various pulsars 
in the Parkes Pulsar Timing Array (PPTA). We then deter¬ 
mine criteria for selecting ‘optimal’ pulsars in section [5] and 
conclude in section [6] 


2 POWER SPECTRUM OF THE PHASE 
RESIDUALS 

2.1 Timing noise 

Let $tn (/) denote the Fourier transform of the autocorre¬ 
lation function of the phase residuals, 5ij>{t), viz. 


/ OO 

dr e 2 ™ fT {S(j>(t)8<j>(t + r)) . 

-OO 


(1) 


If the timing noise is stationary, ( 84>{i)5(j){t + t)} is indepen¬ 
dent of t, as is the mean-square phase residual 


1 r°° 

= df 4? T n(/)■ 


( 2 ) 


In practice, the time spent observing the neutron star, T 0 b s , 
is finite. Hence, one must replace the lower terminal of the 
integral in the right-hand side of © by fobs — 1/Tobs- In 
reality, fitting models to timing data implies PTAs are sensi¬ 
tive to / < fobs [see Coles et al. | (2011 [ | and va n Haasteren fe 
Levin ( 2013| ) for details of timing-model fits in the presence 
of red noise] implying the lower terminal in § depends on 
the PTA data analysis algorithm, with / < f a b s . 

Millisecond pulsar radio timing experiments measure 
$tn(/) oc f~ q at low frequencies, / < lyr -1 , with q ^ 
0 (e.g., Kaspi et al. 1994 Shannon & Cordes 2010| 


Haasteren et al. |201l[ |Shannon et al. 2013). However, the 
observed power law must roll over below some frequency, / c , 
otherwise equation i> implies divergent phase residuals. To 
capture this phenomenologically, we model the spectrum in 
its entirety by 


4>tn (/) = 


Ao 


(i + p/py ' 2 


+ Aw, 


( 3 ) 


which has the observed large-/ behaviour and is even in /. 
In equation (3|, Atn (with units of time) is the dc power 
spectral density, i.e. 4>tn (/ <C fc ) = Atn, which cannot 
be measured directly in existing data sets (Shannon et al. 


20131. In the regime where 4 >tn(/) oc f q , we can express 


the more commonly used root-mean-square-induced pulsar 
timing residuals, or, in terms of Atn and q as 


5.64 




Atn V /2 ( fobs \ 1/2 ( P v 
10~ 10 yr ) ' 1 - 1 ' ( ns - ( ) 


lyr- 1 


1 ms 


where P is the pulsar spin period. For completeness, we 
include a white noise component, Aw, in equation (3|, which 
is observed in all pulsars, dominates for / > 1 yr - , and 
is the only observed noise component in some objects. The 
white component contributes weakly to setting the minimum 
observation time for gravitational wave detection by PTAs, 
the key concern of this article. 

Equation can be compared against predictions of 
phase residuals from the cosmological gravitational wave 
background, $gw(/)- The reciprocal of the frequency where 
the two curves intersect gives the minimum observation 
time, Plain , required before an individual pulsar becomes 
sensitive to a gravitational wave background, 


$TN (T mi n) = $GW . 


( 5 ) 


Equation © provides a quantitative method for deter¬ 
mining when the gravitational wave signal will dominate the 
timing residual power spectrum. We emphasise that this is 
only an indicative threshold for detection; it is not a substi¬ 
tute for a careful signal-to-noise estimate given desired false 
alarm and false dismissal rates. Cross-correlation search al- 


gorithms look simultaneously at a range in / (e.g., Hellings 

& Downs 1983; Jenet et al. 2005; Anholm et al.[2009 

van 

Haasteren et al. 

2009 

1. For example, our definition ( 

5| is 


equivalent to the boundary between the ‘weak signal limit’ 
and the ‘intermediate regime’ as defined in Siemens et al.| 


(20131. While Siemens et al. (20131 calculate a scaling of 


gravitational wave detection significance with time assum¬ 
ing only white timing noise, they also perform simulations 
with red noise assuming q = —3. A future research project 
is therefore to introduce red noise with and without a low- 
frequency turn-over into the analytic calculations of|Siemens| 
et al. (2013). 


It is likely that the near future will see an increasing 
number of PTA pulsars satisfy the condition $gw(/) > 
$tn (/), and that this will occur before a statistically signif¬ 
icant detection is announced. Equation © and the analysis 
presented in this article therefore provide an important in¬ 
put into the time-scale on which this condition will be met 
by individual pulsars, as a prelude to a cross-correlation de¬ 
tection strategy. 


© 0000 RAS, MNRAS 000, 000-000 
































































3 


2.2 Cosmological gravitational wave background 

Pulsar timing arrays are sensitive to gravitational wave 
backgrounds generated by two cosmological source^] binary 
supermassive black holes and vibrations from cosmic strings. 


Despite the above caveats, a power-law model for the 
characteristic strain spectrum from cosmic strings given by 
equation © with —1 < a < 0.8 is a good approximation for 


the PTA frequency band (Maggiore 2000). The predicted 
range for Agw is quoted in Table [l] 


2.2.1 Supermassive binary black holes 

At binary separations where gravitational radiation dom¬ 
inates the orbital dynamics, the SMBH background is 
parametrised as a power law 


M/) = ^4gw 


/ 


yr _ 


( 6 ) 


with a = —2/3 (Phinney 2001). The normalisation coeffi¬ 
cient, Agw, is the subject of intense debate. We utilise the 
most recent predictions by |Sesana (2013b I and Ravi et al. 
(20151, quoted in tablejl] These two predictions assume that 
gravitational wave emission has already circularized the bi¬ 
nary orbits; at binary separations where energy loss to en¬ 
vironments dominates instead, the SMBH wave-strain spec¬ 


trum whitens (Sesana 2013a Ravi et al. 20141. Whitening 


of $gw(/) at low frequencies increases T m i n . 

The one-sided power spectral density of the pulsar 
phase residuals induced by h c (f) is given by 


‘J'gw (/) = 


hc(f) 2 

12tt 2 P 2 / 3 ' 


(7) 


$gw(/) has units of time and can be compared directly with 
$tn (/) as in equation ©• 

2.2.2 Cosmic strings 

Cosmic strings are topological defects that may form in 
phase transitions in the early Universe and produce strong 
bursts of gravitational radiation, which may be detectable 


in PTAs (Damour & Vilenkin 2000 2001 20051. A cos¬ 


mic string-induced stochastic background of gravitational 
waves is characterised by three dimensionless parameters: 
the string tension, Gp, the reconnection probability, p, and 
a parameter, e, related to the size of loops. The best quoted 
limit of Gp < 1.2 x 10 -8 is derived from PTA limits of the 
stochastic gravitational wave background (van Haasteren 
et al.|20lH 20121, although a more stringent constraint (still 


to be computed), is possible with existing data sets [see 
Sanidas et al. (2013) for projected constraints in the near fu¬ 
ture]. Combined observations using the ground-based Laser 
Interferometer Gravitational Wave Observatory (LIGO) and 
Virgo constrain the e - Gp plane to be 7 x 10~ 9 < Gp < 


1.5 x 10 ' and e < 8 x 10 1 (Abbott et al. 


2009 


Aasi 


et al. 2014). Limits on Gp are model dependent; the re¬ 


connection probability is inversely proportional to $gw(/), 
and smaller values of e increase the minimum gravitational 
wave frequency emitted. This can take the maximum of the 
stoch astic background ou t of the sensitivity band for PTAs 
(e.g., Siemens et al.|[2007 Olmez et al. ||20 1*0 i. 


1 Relic gravitational waves from inflation, such as those purport¬ 
edly seen by the BICEP2 experiment (Ad e et al.| [2014), are ex¬ 
pected to be undetectably weak in the pulsar timing band, but 
may be relevant for Advanced LIGO; see |Aasi et al.| ( |2014[ ) and 
references therein. 


3 MODEL INDEPENDENT MINIMUM 

OBSERVATION TIME 

To attain adequate sensitivity to gravitational waves at a 
frequency, /, in the phase residuals of an individual pulsar, 
we must have 4>gw(/) > 4?tn (/) for that pulsar, subject to 
the caveats regarding specific data analysis algorithms ex¬ 
pressed in the text following equation |5|. If equation © 
applies across all relevant frequencies, and <Ftn (/) turns 
over below f c , then 4?gw (/) > $tn(/) is always satisfied 
for some / = as in equation ©. 

In figure [l] we plot 4?tn and 4>gw as functions of /. 
The coloured shaded regions and the region enclosed by 
the black dotted curves in the left-hand plot contain all the 
$gw curves in the parameter range in Table [l] for SMBHs. 
The blue shaded region is the 95% confidence interval from 
Ravi et al. (20151 as described in table [l] The black dot¬ 


ted curves enclose the 68% confidence interval from ISesanal 
(|2013 b Th e shaded grey region is the predicted range from 
Ravi et al. (20141 that includes low-frequency-whitening of 
$gw(/) due to non-circular binaries. In the right-hand plot, 
the green shaded region is the parameter space enclosed by 
the cosmic string predictions from table [T] with a = —1. 
The dotted black curves are specific, representative calcu¬ 
lations of the cosmic string background with p = e — 1 
and Gp = 1.2 x 10 -8 (top curve) and 1.0 x 10~ 10 (bottom 
curve rl The black, red and purple curves in each panel are 
indicative examples of pulsar timing noise as described by 
equation ([3]), with values of f c , Atn, q and Aw given in the 
caption to figure [T| The correspondingly coloured dashed 
curves extrapolate backwards the power-law scaling (equiv¬ 
alently assuming f c —> 0). Finally, the dot labelled ‘PPTA’ 
in both panels marks the lowest limit on the stochastic grav¬ 
itational wave background from Shannon et al. (2013 ). 


Figure [l] illustrates the principal idea of this paper. If 
4>tn is a simple power-law without a low-frequency turn¬ 
over, and for moderate values of q, timing noise masks the 
gravitational wave background down to low frequencies, and 
Tniin is correspondingly long (we quantify this below). A 
turn-over in "Few (/) at some f c is therefore critical for prac¬ 
tical PTA experiments with any millisecond pulsar that ex¬ 
hibits a steep timing noise spectrum with q > 2. The low- 
frequency plateau in 4?gw from elliptical binary SMBHs (the 
grey shaded region in the left-hand panel of figure |T]) makes 
the need for a turn-over in 4?tn (/) even more acute. 

In figure[2]we plot the minimum observation time, T m i n , 
defined by equation ©, as a function of the asymptotic 
(high-/) timing noise spectral index, q, for a hypothetical 
pulsar with P = 10 ms and various values of f c and Atn 
in each panel. The shaded regions and dotted black curves 
delineate the ranges of T , m in(?) for binary SMBHs and cosmic 


2 Calculations used the GWPlotter 
http://homepages.spa.umn.edu/~gwplotter 


website: 
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4 Lasky, Melatos, Ravi & Hobbs 



Source 


A min 
^GW 

A max 

a gw 

a 

SMBHs 

(|Sesana|2013b 

68%) 

3.5 x 10' 16 

1.5 x 10 -15 

-2/3 

SMBHs 

Ravi et ai. 201 

5 95%) 

5.1 x 10 -16 

2.4 x 10“ 15 

-2/3 


Josmic Strings 


1 0 - 16 

io- 15 

-1 < a < -0.8 


Table 1. Theoretical spectral parameters for the gravitational wave background from supermassive black hole binaries (SMBHs) and 
cosmic strings. Three predictions for the SMBH population are presented: the 68% confidence interval from |Sesana |( |2013b| ), and the 95% 
confidence interval from [Ravi et al.] ( |2015| ). The cosmic string models are from Maggiore (2000). 



/ [Hz] / [Hz] 


Figure 1. Power spectral density of timing noise phase residuals, $tNj and gravitational wave phase residuals, $gw (both in s) as 
functions of frequency, / (in Hz). In both plots, the black, red and purple curves are indicative examples of pulsars described by equation 
|3j) with (/ c , Atn, q) = (O.lyr -1 , 10 -10 yr, 4), (0.2yr _1 , 10 -10 yr, 2) and (0.2yr _1 , 10 _8 yr, 2) respectively and with a nominal white 
noise component of = 10 —5 s. The solid curves include the turnover at f c in equation |3]), and the dashed curves extrapolate 
backwards the high-/ scaling oc f~ q (i.e., f c —> 0). The shaded regions and the dotted black curves encompass the regions covered by 
$GW f° r the parameter ranges quoted in Table |7| where the left-hand plot is for SMBHs and the right plot is for cosmic strings. In the 
left plot, the shaded blue region is the 95% confidence interval from [Ravi et al.| ( |2015[ ), the region enclosed by the black dotted curves is 
the 68% confidence interval from [Sesana ( 2013bj) and the shaded grey region is the predicted range from [Ravi et al.] ( |2Q14| ). In the right 
plot, the shaded green region represents cosmic string models with 10“ 16 < Aqw < 10 -15 and a = — 1. The dotted black curves have 
p = e = 1 and Gp = 1.2 x 10 -8 and 1.0 x 10“ 10 for the top and bottom curves respectively. The black dot labelled ‘PPTA’ is the lowest 
published limit on the stochastic background ( [Shannon et al. 2 013 [ l. 


strings, following the same colour scheme as in figure [l] and 
as detailed in the caption of figure [2] The coloured dashed 
curves give the limits on T m i n (q) if 4>tn does not turn-over 
(i.e., f c —» 0). The horizontal dashed black line marks the 
PPTA observing time of 11.3 yr used for the lowest limit on 


the stochastic background published to date (Shannon et al. 


20131. 


To help interpret figures [l] and [2] consider a hypothetical 
pulsar with Atn = 10 -8 yr (i.e., the two left hand panels) 
and q = 2. If the timing noise spectral density turns over at 
f c = 1/5 yr -1 or 1/20 yr -1 , the minimum observation time 
given the most optimistic scenario from Ravi et al. (20151 is 
T min = 20.5 yr or 37.2 yr, respectively. On the other hand, if 
$tn (/) does not turn over, then the dashed blue curves show 
that the pulsar is insensitive to a gravitational wave signal 


until T m in = 70.7 yr. The effect of a plateau in $tn(/) is 
therefore quite striking. Pulsars without a plateau and q> 3 
(depending less sensitively on Atn) are relatively inferior as 
a tool for detecting gravitational waves. 


4 TIMING NOISE FROM SUPERFLUID 
TURBULENCE: A WORKED EXAMPLE 


In section [3j the description of timing noise is model inde¬ 
pendent, in the sense that $tn is parametrised phenomeno¬ 
logically by equation Q, without reference to a specific un¬ 
derlying, physical model. In this section, we repeat the anal¬ 
ysis in section [3] for the timing noise model of |Melatos fc| 
Link (20141 andlMelatos et al. (20151, which attributes the 
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Figure 2. Minimum observation time, T m i n , defined by equation as a function of the high-frequency timing noise spectral index, q , 
with different values of f c and Atn hi each panel. The shaded blue and grey regions represent the ranges of solution space for T m ; n (g) 
for binary SMBHs from Ravi et al. (2015) and Ravi et al.| (2014) respectively, while the dotted black curves encompass the solution space 
predicted by |Sesana| ( |2013b[ l. The shaded green regions are cosmic string predictions. The coloured dashed curves are the corresponding 
limits where the timing noise spectrum does not plateau at / < f c . The horizontal dashed black line is longest the PPTA observing time 
in published data, viz. 11.3 yr ||Shannon et akj2013|. 


fluctuating phase residuals to shear-driven turbulence in the 
interior of the neutron star. We emphasise that we do not 
express any theoretical preference for this model ahead of 
other models in the literature (see section [l]). We focus on 
it here only because (i) it is predictive, (ii) its results can be 
expressed in compact, analytic form and, (iii) the theoretical 
formula for 1 J)tn(/) depends on just three internal neutron 
star parameters, so it is easy to infer constraints on these 
parameters by combining the model with data. 


Consider an idealised neutron star model in which the 
rigid crust is coupled to the charged electron-proton fluid 
which, in turn, couples through mutual friction to the invis- 
cid neutron condensate. The electromagnetic braking torque 
creates a crust-core shear layer that excites turbulence in 


the high-Reynolds number superfluid {Peralta et al. 2005 


2006a|b| Melatos & Peralta 2007 Peralta et al.| 2008|. The 

turbulent condensate reacts back to produce angular mo¬ 
mentum fluctuations in the crust, which are observed as tim¬ 
ing noise ( Greenstein|1970 Melatos & Peralta 20101. In par¬ 
ticular, Melatos et al. (20151 showed that the timing noise 


spectral density can be expressed as 


87 r 1 / 2 r [(q — l)/ 2 ] A 2 r ? (R- 1 ) 

L 


dx x 35 ^ 3 (a ; 4 + 3a : 2 + 9) 


1 + 


4t r 2 / 2 


r/ (i ? -1 ) 2 a: 4 / 3 . 

where V(x) is the Gamma function. Equation | 8 | contains 
three free parameters: the non-condensate fraction of the 
moment of inertia, A = 7 c //o, the decorrelation time-scale, 
77 ( J R“ 1 ) 1 , and q. Here, I c is the moment of inertia of the 
crust plus the rigidly rotating charged fluid plus entrained 
neutrons, Jo is the total moment of inertia, and we define 
r) (-R -1 ) = (27r)~ 1 ^ 2 e 1 ^ 3 7? _2 ^ 3 7, where e is the energy dissi¬ 
pation rate per unit enthalpy (which, in general, is a function 
of the spin-down rate), 7 = Teddy/Tt ur b ^ 1 is the ratio of 
the eddy turnover time-scale to the characteristic time-scale 
over which turbulent structures change (which is longer in 
general due to pinning), and R is the stellar radius. 

The value of the exponent, q, in equation Q depends on 
the form of the superfluid velocity two-point decorrelation 
function. Melatos & Peralta (20101 executed a first attempt 


-q/2 

(8) 
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6 Lasky, Melatos, Ravi & Hobbs 


to calculate the velocity correlation function numerically on 
the basis of Hall-Vinen-Bekarevich-Khalatnikov superfluid 
simulations (Peralta et al. 20081, but it is not well under¬ 


stood for terrestrial turbulence experiments, let alone for a 
neutron star interior, especially when stratification plays a 
role (e.g., Lasky et al. 2013, and references therein). An em¬ 
pirical choice is therefore made that reproduces the asymp¬ 
totic power-law dependence from timing noise data, i.e., 


$tn oc f q as / —» oo [for details see Melatos & Link 


(20141; Melatos et al. ( 2015|)]. We emphasise equation (|8j) 


is not a unique choice, nor can it be inverted uniquely to 


infer the underlying velocity correlation function (Melatos 
et a l. 20151. 


In addition to the power-law scaling at high-frequencies, 
the superfluid turbulence model predicts a plateau at / < 
f c « For time intervals greater than ~ 1 /ri(R~ 1 ), 

turbulent motions throughout the star decohere, implying 
torque fluctuations exerted on the crust become statistically 
independent. By expanding equation | 8 ]l for / <g f c and 
f ^ fc, and evaluating the resultant expression in terms of 
equation CD , we find 

9r(«/2) 


Atn = 


fc = 


16{2n) 67 / 6 V2T[(q - l)/2]r 1 (R- 1 )X 2 


167T 4 + 


120tt 2 


vdr 1 ) 

(27r) 1 / 3 


- 10 


13 
I671" 


45 
+ -tt 


12tt 2 


g-10 9-13 


+ 


9 


167T 4 + 


1207T 2 

13 


45 

+ TT 


9-16 
1/9 


(9) 


( 10 ) 


Equations ( 9J) and (10 1 relate the phenomenological model 
in section [3 to the specific physical model in this section. A 
similar approach applies equally to other models. 

I 11 figure [3] we show four examples of millisecond pul¬ 
sar phase residual power spectra measured by the PPTA 
(Manches ter et al.||2013 |. Overplotted on the data are rea¬ 
sonable ‘by-eye’ fits generated by the superfluid turbulence 
model for 9 = 2, 4 and 6 . The fits are neither unique nor 
optimal (e.g., in a least-squares sense), but they are repre¬ 
sentative. It is outside the scope of this paper to extract 
detailed fits and values for A, r/^Rr 1 ), and 9 for each pul- 
sai[^] We simply note that a broad range of parameters fit the 
phase residuals for any given pulsar. The pulsars shown in 
figure [3] have been chosen as they appear to have moderate 
to high levels of timing noise, cf. other PPTA pulsars. All ex¬ 
hibit a relatively red spectrum. In the context of superfluid 
turbulence, they imply f c > 10 -2 yr _1 , so that the plateau 
is potentially observable in the not-too-distant future]^] 

I 11 figure [4] we plot two further examples of millisecond 


3 The amplitude and spectral index of red-noise in pulsar timing 
residuals are highly covariant, especially when only the lowest few 
frequency bins show evidence for red noise (e.g., |van Haasteren| 
|et al.|2009| |van Haasteren fc Levin|2013| l . Finding best-fit param¬ 
eters for the superfluid turbulence model is therefore a non-trivial 
task that will be the subject of future work. 

4 We note that PSR J1824—2452A resides in a globular cluster 
l |Lyne et al. 1987|, implying most of the timing noise is likely a 
result of motions within that cluster rather than superfluid turbu¬ 
lence. The curves shown in figure [3] therefore represent an upper 
limit on the contribution from superfluid turbulence. 



Figure 3. Phase residual power spectra, < T > (/) for four millisec¬ 
ond pulsars from the Parkes Pulsar Timing Array (thick black 
curves). Overplotted are theoretical curves generated by the su¬ 
perfluid turbulence model for the physical parameters (A, rj, q) 
specified in the legend. The corresponding values of Atn and f c 
for each curve are shown in figure [5] 



Figure 4. Same as for figure [ 3 ] but for PSRs J0437-4715 and 
J1909-3744, two pulsars from the PPTA catalogue with the lowest 
level of timing noise. Theoretical curves for the superfluid turbu¬ 
lence model require smaller values of f c //(II 1 ) than in figure 
[ 3 ] to remain consistent with the data. 


pulsar phase residuals. These objects exhibit the lowest level 
of timing noise in the PPTA sample. For the superfluid tur¬ 
bulence model to remain consistent with these data, the ob¬ 
jects must have long decorrelation time-scales, i.e., f c < 10 2 
yr _1 . The data show the white noise component, Aw, and 
the turbulence-driven red-component sits below Aw- Under 
these circumstances, the turnover in <I>tn(/) occurs too low 
in frequency to be observed, and the main factor limiting 
PTA detection is Aw- 


© 0000 RAS, MNRAS 000, 000-000 























































































7 


5 OPTIMAL PULSARS 

What pulsars are best placed to detect a gravitational wave 
background, given the longest time one is prepared to wait? 
In figure[5]we plot 1 / f c against Atn, for different values of q 
and Tinm in each panel. The left-hand vertical axis displays 
the results for the model-independent form of $tn in equa¬ 
tion The right-hand vertical axis registers the decor¬ 
relation time 1 / 77 (T 2 1 ), in the superfluid turbulence model 
in section [4] The dashed grey curves are curves of constant 
A. Overplotted are the superfluid turbulence model ‘fits’ to 
the PPTA pulsar data in figure [3] where the open circles, 
closed circles, open squares and closed squares are PSRs 
J1024—0719, J1643—1224, J1824-2452A and J1939-2134 
respectively. 

Figure [5] allows us to ask whether, for example, 20 
yr of timing a specific pulsar will allow for sensitivity to 
the most optimistic SMBH gravitational wave strain of 
Agw = 2.4 x 10 -15 . I 11 the middle set of panels, the latter 
strain limit appears as the right-most boundary of the blue 
shaded region. A pulsar with timing noise below this curve 
is sensitive to a gravitational wave signal in T 0 b a ^ 20 yr. 
Sensitivity depends on q as illustrated in the three different 
panels running vertically. It also depends on f c . For exam¬ 
ple, a hypothetical pulsar with q = 4 and Atn ~ 10 -11 
yr is only sensitive to a gravitational wave background for 
1/fc < 12 yr. This is an interesting constraint: a pulsar in 
a PTA that tolerates T m i n ^ 20 yr is sensitive to a gravita¬ 
tional wave background if 4>tn(/) exhibits a plateau after 
< 12 yr of timing. 

The superfluid turbulence model fits from figure [3] give 
an indication as to the usefulness of individual pulsars from 
the PPTA dataset. For example, consider PSR J1939—2134 
(closed squares). If one again tolerates T m i n ^ 20 yr, the fits 
imply a pulsar is sensitive to a conservative prediction for 
the gravitational wave background for q < 4, although for 
q ~ 4 this requires the timing noise spectrum to plateau after 
approximately 15 yr of timing. We emphasise again that the 
model fits should only be taken as indicative; careful and 
detailed analysis is required to extract the true timing noise 
signal parameters from the data. 


6 CONCLUSION 


Pulsar Timing Array limits on the cosmological gravita¬ 
tional wave background are continually dropping to the 
point where they usefully constrain galaxy formation mod¬ 
els (Shannon et al. 20131. Positive detections, on the other 
hand, require a cross-correlation algorithm to simultane¬ 
ously analyse timing residuals from multiple pulsars. Such a 
detection will likely occur when the gravitational wave back¬ 
ground is the largest component in the unmodclled portion 


of many individual pulsar’s timing residuals (Siemens et al. 


20131. If the timing noise spectrum is steeper asymptoti¬ 
cally (at high /) than the gravitational wave spectrum, this 
is only possible if the timing noise spectrum flattens below 
some frequency, f c . In this paper, we calculate the minimum 
observation time required, given / c , before the gravitational 
wave background rises above the timing noise plateau in any 
specific pulsar. We calculate this minimum observation time 
both in a model-independent way, and for timing noise aris¬ 
ing from superfluid turbulence. The latter model is selected 


not because it is necessarily preferred physically, but be¬ 
cause it is simple, predictive and analytically tractable and 
therefore provides a test-bed for repeating the calculation 
with other physical models in the future. 

Our results rely on the timing noise spectrum whiten¬ 
ing below some threshold frequency, / c . This provides an 
observational diagnostic that can be used to infer the effec¬ 
tiveness of an individual pulsar in a PTA. If, upon observing 
a pulsar for some T > 1 /f c , one finds that <1>tn (/) has not 
whitened below / c , that pulsar’s capacity for assisting use¬ 
fully in the detection of a gravitational wave background is 
severely diminished. The f c for a given pulsar is a function 
of the rotational parameters of the pulsar, and the gravita¬ 
tional wave amplitude and spectral index. Therefore, using 
the prescription outlined in this paper, one can predict f c for 
a given pulsar and a given gravitational wave background. 

In reality, measuring f c in a single pulsar is difficult. 
Firstly, the noise in a given pulsar timing power spectrum is 
large, and secondly, the power in the lowest-frequency bin is 
generally dominated by the fact that a quadratic polynomial 


is fit to the timing residuals [see van Haasteren & Levin 
( |2013[ )]. These two effects potentially mimic a low-frequency 
turn-over, implying multiple low-frequency bins are required 
to confirm the existence of a low-frequency cut-off. 

Many data analysis algorithms simultaneously fit the 
timing model and the unknown noise contributions for any 
individual pulsar. In this sense, one can include a low- 
frequency plateau into gravitational-wave detection algo¬ 
rithms, e.g., by way of a Bayesian prior on the form of the 
power spectral density. Physically motivated models for tim¬ 
ing noise, such as the superfluid turbulence model discussed 
herein, could be used to guide such priors. 
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